/************************************************************************
*Description: Program for parallel execution of matrix multiplication	*
*Authors:	  Amber Dhavale, Chitrangda Potdar, Ronald Ngatuni			*
*Version:	  1.0														*
************************************************************************/

//Inclusion of header files
#include <iostream>
#include <cstdlib>
#include <time.h>
#include <omp.h>

using namespace std;

/********************************************
* Name: main()								*
* Description: Program execution starts here*
********************************************/
int main()
{
	double startTime,endTime = 0.00;
	startTime = omp_get_wtime();

	//Declaration of variables
	int const rowA = 5000;
	int const columnA = 2000;
	int const rowB = 2000;
	int const columnB = 2500;
	int const rowC = 5000;
	int const columnC = 2500;
	int num_of_threads = 8;

	//Declaration of matrix variables (pointers for dynamic allocation)
	double **matrixA,**matrixB,**matrixC;

	//Mutex variable
	omp_lock_t mutex;

	//Allocating memory for Matrix A
	matrixA = new double*[rowA];
	for(int rowIndex = 0; rowIndex < rowA; rowIndex++)
		*(matrixA + rowIndex) = new double[columnA];

	//Initializing Matrix A
	srand(time(NULL));
	for(int rowIndex = 0; rowIndex < rowA; rowIndex++)
		for(int columnIndex = 0; columnIndex < columnA; columnIndex++)            
			matrixA[rowIndex][columnIndex] = (rand() % 99 + 1)/10.0;

	//Allocating memory for Matrix B
	matrixB = new double*[rowB];
	for(int rowIndex = 0; rowIndex < rowB; rowIndex++)
		*(matrixB + rowIndex) = new double[columnB];

	//Initializing Matrix B
	for(int rowIndex = 0; rowIndex < rowB; rowIndex++)
		for(int columnIndex = 0; columnIndex < columnB; columnIndex++)
			matrixB[rowIndex][columnIndex] = (rand() % 99 + 1)/10.0;     

	//Allocating memory for Matrix C
	matrixC = new double*[rowC];
	for(int rowIndex = 0; rowIndex < rowC; rowIndex++)
		*(matrixC + rowIndex) = new double[columnC];

	//Initializing Matrix C
	for(int rowIndex = 0; rowIndex < rowC; rowIndex++)
		for(int columnIndex = 0; columnIndex < columnC; columnIndex++)             
			matrixC[rowIndex][columnIndex] = 0;

	omp_init_lock(&mutex);
	int rowIndex = 0; 
	int columnIndex = -1;
	int numberOfElementsInPool = rowC * columnC;
	int elementDoneComputing = 0;
	int i;
	int j;
	omp_set_num_threads(num_of_threads);

	#pragma omp parallel shared(rowIndex, columnIndex) private(i,j)
	{
		while(elementDoneComputing < numberOfElementsInPool)
		{
			omp_set_lock(&mutex);
			if(columnIndex < (columnC - 1))
			{
				columnIndex++;
				i = rowIndex;
				j = columnIndex;				
			}
			else
			{					
					if(rowIndex < (rowC - 1))
						rowIndex ++;
					else
					{
						//cout << "I am thread #" << omp_get_thread_num() << " entered more than once so exiting" << endl;
						omp_unset_lock(&mutex);
						break;
					}
					columnIndex = 0;
					j = columnIndex;
					i = rowIndex;
			}
			//cout << "I am thread #" << omp_get_thread_num() << " working on element C[" << i <<"][" << j << "]" << endl;
			omp_unset_lock(&mutex);
			
			for(int k = 0; k < columnA; k++)
				matrixC[i][j] += matrixA[i][k] * matrixB[k][j];

			elementDoneComputing++;
		}
	}
	omp_destroy_lock(&mutex);
	endTime = omp_get_wtime();

	//Display the time required
	cout << "Total Time required = " << (endTime - startTime) << " seconds" << endl;

	/*for(int rowIndex = 0; rowIndex < rowC; rowIndex++)
	{	
		for(int columnIndex = 0; columnIndex < columnC; columnIndex++)             
		{
			cout << matrixC[rowIndex][columnIndex] << " ";
		}
		cout << endl;
	}*/
	//system("Pause");
}
